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Rapid evaluation of the periodic Green's function in d dimensions 
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Abstract 

A method is given to obtain the Green's function for the Poisson equation in any arbitrary 
integer dimension under periodic boundary conditions. We obtain recursion relations which relate 
the solution in d-dimensional space to that in (d — l)-dimensional space. Near the origin, the 



Green's function is shown to split in two parts, one is the essential Coulomb singularity and the 
other part is regular. We are thus able to give representations of the Coulomb sum in higher 
dimensions without taking recourse to any integral representations. The expressions converge 
exponentially fast in all part of the simulation cell. Works of several authors are shown to be 
', special cases of this more general method. 
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I. INTRODUCTION 



The Poisson equation is probably one of the most useful equations in physics. In a two- 
dimensional (2D) space, the periodic solution of this equation corresponds to the solution of 
particles interacting with the logarithmic interaction, and it has applications in simulations of 
2D pancake vortices in high-temperature superconductors^. In 3D, periodic solutions to the 
Poisson equation are used in electromagnetism. Here, the solution of the Poisson equation 
corresponds to a number of charges interacting with the Coulomb potential. This 3D periodic 
solution is routinely used in most simulations involving charged particles. Recently, the 
periodic solution of the Poisson equation in higher dimensions has found use in the string 
theory. 

In ID and 2D, the Green's function for the Poisson equation for a charge neutral box may 
be obtained in a closed form. In 3D, one can obtain rapidly converging series representations 
using well known method by Ewald^. The other two approaches for the 3D case were given 
Leknep^ and SperfA However, in higher dimensions, one can either use the Ewald method 
which has its drawbacks, or use the Jacobi theta function identities^. In general, there is no 
efficient way to calculate the Green's function in a general (^-dimensional space with d > 3. 

In this paper, we give an exponentially fast converging series representation for the 
Green's function of the Poisson equation in any positive integer dimension. This work 
will generalize the methods employed for 2D and 3D case^, and will tie together the differ- 
ent approaches taken by Lekner 3 and Sperb* for the especial case of d = 3. The outline 
of the paper is as follows. In Sec. I we derive expressions giving the Coulomb sum in the 
ci-dimensional space. In Sec. II we derive recursive relations using the result of the previous 
section. In Sec. Ill we discuss the results. 

II. GREEN'S FUNCTION IN d DIMENSIONS 

For simplicity, we consider the case of a unit charge situated within a cubic box in d 
dimensions. The sides of the box are all assumed to be of unit length. From here onwards, 
we will refer to the box as simulation cell. The basic simulation cell repeats itself in all d 
dimensions. We also assume a charge neutral system. The unit charge interacts with other 
identical unit charges (for the case of different charges q± and q<i one just gets an extra factor 
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of 9152) situated at the vertices of the periodic structure. The periodic Green's function in 
d dimensions satisfies the Poisson equation, 

V 2 d G(r) = -C d J2Kr+l), (2.1) 
1 

where is the Laplacian operator in d dimensions, I denotes a d- dimensional vector, whose 
components are integers ranging over —00 to +00, and C d is a dimension- dependent factor. 
The value of C d for various dimensions is 

2 d=l 
( 'd = < 2tt d = 2, 

4n u+1 /T (u) d>2. 

Here, T(u) stands for the Gamma function, and v = (d — 2) /2. We note that with this 
choice of C d in Eq. (J2.1j) . the G stands for the Coulomb type summation in d dimensions. 
Thus, G corresponds to a sum of type — \r\ in ID, a logarithmic sum, -In | t* | , in 2D and 
a sum of type |r| - ^ -2 ^ for a <i-dimensional space with d > 2 . The solution of Eq. ()2.1|) 
diverges, which is a simple consequence of the fact that the interaction energy of a charge 
with another charge and all its periodic images is infinite. To obtain a meaningful value of 
G we will have to modify Eq. (|2.1j) as follows 6 : 

V 2 G d (r) = -C d 8(r+l) + (2.2) 

The second term in eq. (|2.2|) amounts to the presence of a uniform background charge. Thus, 
for every charge, q, one may imagine a uniform distribution of charge, such that the total 
charge per basic simulation cell adds up to —q. For a charge neutral periodic system, 
imposing these kind of background uniform charge distributions does not matter since the 
total uniform background charge adds up to zero. However, now a unit charge located within 
the basic simulation cell at position {xi} not only interacts with a second charge located 
at the origin and its periodic images, but also interacts with the neutralizing background 
charge of the second particle. This particular way of introducing the artificial neutralizing 
background charge leads to only the intrinsic part 3 of the potential energy. We note that 
once the Green's function is obtained, the solution of the equation 

V 2 V d = -C dP {r) 



under periodic boundary conditions could be simply obtained from 

V d = f G d (r-r')p(r')dr', 

J cell 

where p is periodic and the simulation cell is overall charge neutral. The rapid evaluation 
of the G d is discussed in the next section. 

The solution of eq. (|2.2jl can be written easily in the Fourier space as^ : 



Gd(Xi, X2, 




i2n(m 1 x 1 +m2X 2 A \-m d x d ) ^2 



(2.3) 



\ [m\ + m \ -| \~m 2 d + P 2 /4ii 2 } (3 2 

where (3 is an infinitesimal parameter which tends to zero. Here , the set {rrii td } denotes 
a set of d integers {mi,m2, ...,m d }. Each one of these integers mi runs over — oo to +00. 
Also, Xi,x 2 , ■ ■■,x d denote the components of vector r d in d- dimensions. Due to the periodic 
boundary conditions, it is sufficient to treat the case where each Xi satisfies —0.5 < x < 0.5. 
The complete expression for the potential has a term arising from the surface contribution. 
For the 2D case this term turns out to be zero, but for 3D one obtains a contribution from a 
dipole terrrA At this point, we would recast the Eq. ()2.3|) in an alternative form. For that, 
we use the fact that the solution of 



(V 2 -/? 2 ) Q (r) = -6(r) 



in (i-dimensional space is given by 



1 P"K v (l3\r\) 



Qo(M;/?) = (2 ^ +1 r „ . 
Thus, the solution of 

(V 2 -e 2 ) Q d (r;(3) = -C d J2S(r+l) 

1 

in (i-dimensional space will be given by 



Qd(r;0 



C d 
(2nf 



(2.4) 



(2.5) 



(2.6) 



(2.7) 



{ m l,d} 



l,d 



where 



1/2 



rid 



E 

i=l 



771 i Xo 



(2.8) 
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On the other hand, the solution of Eq. (|2.4j) can be written down in the Fourier space easily 
as 



Q d {r;(3) 



Ca 



i2ir(mixi+m2X2-{ \-m d x d ) 



(2tt) 2 frf {m\ + m\ + • • • + m 2 + /? 2 /4tt 2 } ' 
Using Eqs.(j23J) and (|2.4jl we see that one can write 

/ 

— V 



(2.9) 



G d (x 1 ,x 2 , ■ ■ ■ , x d ) = C d lim 



(27T) 



{ m l,d} 





\ 


K v {(3r ljd ) 

'l,d 


1 







(2.10) 



A yet another alternative form of G d can be obtained as follows. We can perform one of the 
d sums in Eq. fl2.Hj) analytically using the formula 8 

exp (2iTimx) 7r cosh [777 (1 — 2 \x\)] 



E 



Thus, we obtain 



m 2 + 7 2 7 sinh (717) 
/ 



(2.11) 



Gd{x\i x 2 , ■ ■ ■ , x d ) 



., lim 



E 



^ cosh 


( x - 2 




^{m 2 , d } sinh 







x exp 



2ix% mjXj 

i=2 



1 

J 2 I ' 



(2.12) 



where 7{ m2 d } is defined as 



^4= X>* +/52 



1/2 



(2.13) 



,i=2 



For the purpose of taking the limit /3 — *> 0, the sum in the first part of Eq. f|2.12|) is broken 

as 

2^ = /J + (Term with m 2 = 0, m 3=0 ..., m d = 0) , (2-14) 

{"12, d) {™2,d} 

where a prime over the summation sign indicates that the term corresponding to all m ; 
being zero is to be excluded from the summation. This leads to the following representation 
for G d : 



n cosh 


P{m^}0-- 2 \ X l\)_ 


^{»™2,d} 


Sinh (^{ m2 , d }) 



x exp 



2iri m i x i ) + H d 



[2.15) 
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where we have taken the limit (3 — > 0, i.e. we have substituted (3 = in the first part, and 
H d is given by 



H d = Cd n lim 



2tt 2 cosh [(1/2 - \xi\)0\ 4vr 5 



(2tt) 2 ^oV/3 sinh (0/2) P 
C d ^ (l-6|xi|+6x?). 



(2.16) 



To avoid the bad convergence towards x\ — > 0, we further modify the summation in the first 
part of Eq. (|2.15J) by using the following trigonometric identity 



cosh (a — b) cosh(a) . . 

= exp(— b) , - ^ + exp(— a) 



sinh (6) 



sinh (b) 



Thus, can be written as 



G d = H d + J d + M d , 
where H d is defined in Eq. (I2.1fij) . J d is given by 



(2.17) 



(2.18) 



7T 



exp 



x 



{m 2 , d } 1 i m ^} 

cosh [vr7{ m2 d } (l-2|xx|) 



sinh 



exp | 2ni m^Xj 

i=2 



(2.19) 



and 



(27T) 



7T 



exp 



-2 IxJ 7T 



exp 



2 «E 



[2.20) 



i=2 



{m 2 , d } 1 i m ^} 

It is easy to see that Eq. ()2.19|) does not have any convergence problem as attends to zero. 
Thus, the whole problem has reduced to evaluating the M d term efficiently This will be 
done in the next section. 



III. RECURSIVE FORMULAS 

In this section we obtain recursive formulas for G d in two different ways, starting with 
the expressions in Eq. ()2.10|) and (|2.18j) respectively. The first method, with Eq. (|2.10j) as 
the starting point, will contain Lekner's results for d = 3 as a special case, while the second 
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method will contain Sperb's result in 3D as a special case. With the help of Eqs. (j2.10J) and 
()2.3j) we can write 



G d (x 1 ,x 2 , :.,x d 



a, 



(2vr) P-*o 

— ^-s lim 

(2tt) P-+o 



lim Q d {x 1 , x 2 , .., x d ; (5) - — 



( 



exp (2-Ki m i X i^ ^ 

exp (2nim 1 x) x ^ d ^ — 2 



y{^i,d} 

Using the definition of G d , Eq. ()3.1|) be written as 

C d 



jr (3-D 



G d (xi,x 2 , ...,x d 



lim 



J] exp (2 



7T«miXi / 



xQd_i (x 2 , ..x d ; \J [3 2 + (27rmi) 2> ) - 



(3.2) 



We separate out the term corresponding to ra\ = in Eq. ()3.2|) so that the limit corre- 
sponding to (3 can be taken. Thus, we write Eq. (j3.1|) as 

c 00 

G d (xi, x 2 , x d ) = 2— — X~] cos (27rmixi) (x 2 , x 3 , ...,x d ; 2irm 1 ) 



mi=l 



+ 77-— lim 



Qd-i (^2, £3, x d ; j3) - 



2—^- ^ cos(2vrmiXi)Q d _i (x2,X3,...,x d ;27rmi) 



mi=l 



+ 7^ — Grf_i(a;2, .., x d ), 



(3.3) 



where we have taken the limit (3 — > in the first term. The Eq. (|3.3|) is one of the most 
important result of this paper. This relates a d- dimensional sum to a (d— 1) -dimensional 
sum. This is a recursive relation. If one is able to obtain the Green function for the (d — 1)- 
dimensional space, one can obtain the Green's function for the d-dimensional space. The 
first term in Eq. (J3.3)) can be modified in the following way. We can use a form of G d _i 
similar to the one used in Eq. (|2.7|) to obtain 



G d (xi,x 2 , x d ) — 2 



C d 



(2tt 



1/+1/2 



2_j c °s (2nmiXi 



m\=l 



{ m 2,d} 



x (27rmi 



v-i/2 K v-i/2 (27rmir 2 , d ) 
' u-1/2 



C d 
C d -\ 



G d -i (x 2 , ..,x d . , 



(3.4) 
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where {m 2 ^} denotes a sum over sets {m 2 m 3 , ...,ra<z} and r 2t d is defined like Eq. (|2.8|) 

1/2 



?~2,d 



(mi - XiY 



i=2 



(3.5) 



Let us now consider three different cases corresponding to d — 1, d — 2 and d > 2. For = f 
we can evaluate Gd=i in a closed form: 



Gi 



(2tt) z /»-o 
' lim 



7711 



exp (27iimiXi) 1 
£ 2 + m? ]P 



7r cosh [7T/3 (1 — 2 \xi\ 



(2-k) P-+o\/3 sinh [tt/5] 
Ci^(l-6|xi| + 6a^). 



1 

J 2 



Also, the self-energy for this case may be obtained as 

Gf { = lim Gi + = 



zi^O 



Ci 
12' 



For d — 2 case, we obtain using Eq. ()2.1()j) : 

„ oo +oo 

C-2 



G 2 (xi,x 2 ) = 2 ^ cos (2-KmiXi) (2nmi 

(2tt) ' — ' — ' 

V / 77ii=l m2=— oo 

K_i/ 2 (27rm! |x 2 + m 2 |) , C 2 



-1/2 



\x 2 + m 2 \ 



-1/2 



—G 1 (x 2 ) . 



Now, using the relation 8 



K -v* ( r ) = \/^ ex P 



we can write 



C 2 (a;i,a; 2 ) 



C 2 cos (27rm 1 x 1/ 

2tt 



E E 

7712= — OO 7711=1 



mi 



C 2 



x exp (-27rmi |x 2 + m 2 |) + ttGi (x 2 

The sum over mi can be easily carried out using the identity- 

x cos(27rmiXi) 
L(xi,x 2 ) = 2^ exp (-27rmi |x 2 |) 



7711 =1 



(3.6) 



(3.7) 



(3.8) 



(3.9) 



In (1 — 2 exp [— 27nr 2 ] cos [27rxi] + exp [— Anx 2 }) . 



(3.10) 
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Thus, Go can be written as 



Q +00 

G 2 (xi, x 2 ) = 77^ L(xi, \x 2 + m 2 |) + L(xi, \x 2 - m 2 \) 

m 2 =l 

+ L(:r 1 ,x 2 ) + - 2 -G' 1 (x 2 ). (3.11) 



It is also trivial to derive 



G self = 2 g £ L(Q) |m2 | } _ + ^ (3 12) 

m2=l 



Now we consider the case for d > 2. We can obtain from Eq. ()3.4j) . It is seen that for 
large arguments the modified Bessel functions decay as 

Mr)~y|Jexp(-r). (3.13) 

As a result, the first term in Eq. ()3.4j) decays exponentially. However, one may run into 
problem if r 2 ,d is very small. In such a case the terms corresponding to {m 2 ,d} all being zero 
form a very slowly converging series over m\. This problem of slow convergence when r 2 , d 
is small can be handled in the following recursive manner. We separate out the particular 
terms corresponding to {m 2 ^} all being zero, and define 

c 00 

E d {xi , x 2 , Xd) — 2- +1 y 2 2_j cos(27rmiXi) 
(2?r) "* =1 



x (27rmi) 



2,-1/2 Ky-i/2 (27rmir 2;d ) 



r 2,d 



I/-1/2 

/ 

c d 



, cos(27rmix/ 
(2tt)^ +1/2 ^ ^ 



x (27rm.i) 



1 { m 2,d} 

2,-1/2 Kv-1/2 (27rmir 2id ) 



„f-l/2 



c 

(2tt) j 



x ( 27mii r 1/2 ^-V^mxr) 



(3.14) 



Now, we show how to handle the evaluation of Ed corresponding to d > 3. The case for 
d = 3 will be almost the same. Using the relation^ (which by the way can be derived from 



Eq. (EH)) 



E 



.2i/ 



X 



VK V (2vr/r) cos (2tt/x) I v > 0, (3. 15) 



we can write 



E d (xi,x 2 , ...,x d ) = 2- — cos(27rmiXi) 



(27T) 

(27rmi) 

oo 

+ E 



= 1 {™2,d} 



Jy _ 1/2 ^-1/2 (27rmir 2 ,, 



r "-l/2 
2,d 



1 



v^F r(z/-l/2) 



A __^[(s + *) 2 + r»] r( 
Also, the sum over in Eq. (J3.16|) can be written as 



V) r 



2i/-l 



(3.16) 



E 

fc=— oo 



1 



[{x + kf + r 2 ] V (x 2 + 



i 



+ 



+ fc) 2 + r 2 ]" [(ar-jfef + r 2 ]" 



+ £ ( !'\ 2l [((2l + 2v,N + x) + C (2/ + 2z/, AT - 



£)], (3.17) 



where A/" is an arbitrary integer^ such that N > r + \x\. Using Eqs. (|3.4j) , (|3.16|) and (|3.17|) 
we can now write 



G d (xi,x 2 , ..x d ) 



1 



C d 



cos (2nmiXi) 



x (27rmx) 



„_ 1/2 /C-1/2 (27rmir 2 ,, 

I/-1/2 



2,d 



JV-1 



+ E 



+ 



+E 

i=i 
C d 



" x \[{x! + kf + rlY - kf + r]]\ 

rf [C (21 + 2u, N + x x ) + C (2/ + 2z/, iV - Zl )] 



+ 



,2;/ 



(3.18) 
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Note that if d = 3 then instead of Eq. ()3.17|) we should use 

oo 

4 ^2 K o (2-Kmx [x\ + xl) 1/2 ^j cos (27rmi 



mi=l 



l 1 X lj 



(xi + xl) l/2 \ I i 




2 <{ 7 + In ( v 2 J/ ) }>+ , 2 = + S(si,s 2 ,2: 3 ), (3.19) 



where 

S (x x ,x 2 ,x z ) = ^2 I + / 

n=i \ y zjj + ac| + (n + xi) 2 y + z§ + (n - 2i) 2 # 

+ 1 - 2 7 - [V(iV + x x ) + </>(iV - a*)] 

+ E ( _1 / 2 ) ( x 2 + [C (2/ + 1, N + x) + C (21 + l,N-x)\. (3.20) 

Thus, for the 3D case one would make the following two changes in the expression given in 
Eq. f 13.18)) . First, there would be an extra term containing — 27 — [if)(N + x\) + ip(N — x\)\ 
on the right hand side, and second the last term in Eq. (|3.18j) would be changed to 

-^[G rf _ 1 (x 2 ,x 3 ) + lnrJ. (3.21) 

Eq. ()3.18j) provides us with a general algorithm to calculate Gd efficiently in any dimensions. 
For an example, if we had started out with d = 10, we can obtain G 10 — r$ 8 by calculating 
Gg— 7vT 7 . Continuing in this fashion we will come down to calculating G 2 +ln r 2 . Now, this last 
part G 2 + lnr 2 has been obtained by several authors . In fact, it can be obtained in a closed 
form 5 . Thus, we have been able to calculate Giq — r$ 8 from which can obtain G10 by taking 
the radial part r$ s on the other side. Other forms of G 2 are given by Gr jb nbech- Jensen^ 
and Tyagi^. For the sake of completion we write down the result for G 2 . 



„ , . 1 v-^ n exp (— 7r Iml) cosh [27rmxi] . . 
G 2 {x 1 ,X2) = — y -, — r -WI — f; cos(27rmx 2 ) 

271 ' ' \m \ am h I ir \m \ \ 



in 



71 


exp (— 71 


m\) cosh [27rmxi] 


\m 


sinh (tt 


m\) 



— — In [cosh (2ttx\) — cos (27rx 2 )] 

+ I(l + 8 ^-!5M. (3.22) 



In the closed form G 2 is written as^ 

(x 2 In 2 1 
Z D7T 2,71 

11 



■&X [tt (xi + ix 2 ) , exp (-71-)] 



^[0,exp(-7r)] 



1/3 



(3.23) 



where i?i represents the Jacobi theta function of the first kind. Also, the self-energy for the 
2D case can be obtained from Eq. ()3.22|) 



G 



self 



oo 

7T ^ 

m=l 



7i 


exp (—7i 


\m ) 


\m 


sinh (j[ 


m ) 



- In (2tt) + 



7T 

6' 



(3.24) 



or it can be obtained from Eq. ()3.23|) : 



lf In 2 2 
G s f = lnvr--ln 

2 3 3 



tf'i(0,g) 



(3.25) 



All three forms Eq. ()3.12j) , Eq. ()3.24j) and Eq. 1)3.25)1 are equivalent and give numerically 
the same value for the self-energy. Similarly Eqs. (j3.ll)) . (j3.22j) and ()3.23)) show perfect 
agreement. 

Now, we give another alternative approach. This time we start with Eq. ()2.18)) . where 
H d , J d and M d are defined in Eqs. ()2.16|) . ()2.19|) and ()2.2(Jj) . H d and J d do not have any 
convergence problem in the region of interest. We show how to handle M d . A recursion 
formula similar to Eq. ()3.4|) can be established for M d . It is easy to see just by inspection 
that M d obeys the following recursion formula: 



\l, , -J-M d „ x + 2- — -4tt^ cos ( 27rm 2£ 2 ) (2vrm 2 )^ 1/2 

m,2=l 



(2n) u+1/2 



K v -!/2 (l-Km 2 \Jx\ + (m 3 - x 3 ) 2 + .. + (m d - x d ) 2 \ 
£ r ,,-1/2 —. < :1 ^> 



{ m 3,d} 



\Jxl + {m 3 - x 3 ) 2 + .. + (m d - x d y 



where M d _ ly analogues to Eq. (|2.2(Jj) . stands for 
C d . 



(2tt 



■i E ^- 



■ exp 



-2 xi 



7r '^{m3,i} 



exp 



i=3 



{-3,4 1 i m ^} 

In the final step, we break the sum in the second part of Eq. ()3.26|) as follows 

E = E+ E ■ 

{m'i,d} {™3,d} m 3 =0,m 4 =0,.. 



(3.27) 



(3.28) 
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The term corresponding to m 3 = 0, = 0... gives rise to a term F d in Eq. (|3.26|) : 

n x 

7-. o 



(27T) 



f+1/2 



cos (27rm 2 X2) (27rm 2 ^ 



m,2=l 



X 



^-1/2 (2itm 2 ^/x\ + xjj +T7Tx|j 



A/xf+xl + + X 2 



■1 w-1/2 



1 



a/xi + x| + xj + .. + x 2 d C d -i yj~x\ + x\ + .. + x\ 
JV-1 



77 + 



J V [(x 2 + fc) + r 2 ] [(x 2 - fc) + r 2 ] 



— v 



r M [C (2/ + 2v, N + x 2 ) + ((21 + 2u, N - x 2 )\ . 



(3.29) 



Thus, we finally obtain the following recursion relationship for M d . 



M d - 



1 



C, 



(x\ + r 2 ) 



2^+1/2 



Cd-i 



+ 2 , n J+i/2 cos ( 27rm 2^2) (27rm 2 )^ 1/2 
(2vr) 



7712=1 



{ m 3,d} 



/ ^-1/2 {2Txm 2 \Jxl + (m 3 - x 3 ) 2 + ... + (m d - x d ) 2 ^j 



^/x 2 + (m 3 - x 3 ) 2 + ... + (m d - x d y 



v-l/2 



JV-1 



77 + 



^ \ [(x 2 + fc) 2 + r 2 ]' / [(x 2 -fc) 2 + r 2 ] 1/ 



+ £ 



rf [C (2/ + 2z/, iV + x 2 ) + C (2/ + 2i/, N - x 2 )} , (3.30) 



where 



+ x 3 -|- ... -|- x^. 



(3.31) 



For d = 3 case, once again, we will have to make two modifications in Eq. (J3.3Uj) . With 
this approach we have obtained Eq. ()3.3()|1 . which is analogues to Eq. ()3.18j) . However, 
the analysis has become a little bit tedious. The advantage of the second method is that it 
reduces the computation time, as there is one less summation. The second advantage it can 
be written down in a product decomposition form. For example, how such a product decom- 
position form may be written, one may consult Sperb, where a special case corresponding 
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to d = 3 is considered. In general, the procedure of dimensional reduction is to be continued 
until we have Mi on the left hand side. It is clear that Mi = 0. Let us again consider three 
special cases. For d — 1 one only has H d=1 and thus G\ — Hi. For d = 2 one obtains 

x " sh H r (l-2 M)] 
smh(7T7 m2 ) 

and M 2 from Eq. ()2.2U|) and (j3.10J) turns out to be just L (x±, x 2 ) 

M 2 = L{x 1 ,x 2 ). (3.33) 

Combing H 2 , J 2 and M 2 we obtain the form of G 2 given in Eq. ()3.22|) . Considering finally 
the case for d > 2 case, we can obtain Gd again from Eq. (|2.18jl . Now and Hd are 
convergent and M^ can be obtained using the recursive relation Eq. (J3.3U)) . For example: 



x 



^i^o ^2irm 2 ^ x\ + (m 3 - x 3 ) 2 ^ 



W3 

iV-1 



+ £ V [(a* + kf + ri] 1/2 + [(x 2 - fc) 2 + ri] 1/2 . 
- 2 7 - [ip(N + x 2 ) + ^(AT - x 2 )} 

+ J2 ( 7 ) r ? K (2/ + 2i/, iV + x 2 ) + C (2/ + 2zv, iV - x 2 )] , (3.34) 
2=1 ^ ' 

where M 2 has already been evaluated above. We see that in all the case, expression could be 
written in a form that the essential Coulomb singularity as the two charges approach each 
other has been removed. 



IV. CONCLUSIONS 



Using the limiting behavior of the modified Bessel functions, we showed how conditionally 
convergent Coulomb sums may be handled in an elegant way. We gave two representations of 
the Green's function for the Poisson equation in any integer dimensional space. A recursive 
method was derived that can be applied for wholly periodic cases, as well as for those cases 
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where one may have open boundary conditions along one of the directions. The method may 
be extended to cover the case where any number of directions may be open. The formulas 
obtained show rapid convergence in all part of the simulation cell. This method is general 
enough that it can be easily generalized for a higher dimensional "triclinic" cell. A particular 
case of the application of this method for a triclinic cell can be seen in a recent paper—. We 
have shown that the present work generalizes the work of several authors on periodic and 
partial periodic systems 3 -! 4 *^. To our knowledge, this treatment is the first of its kind ever 
taken in a dimension higher than d = 3. 
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